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We calculate the quadratic fluctuations of net baryon number, electric charge and strangeness 
as well as correlations among these conserved charges in (2+l)-flavor lattice QCD at zero chemical 
potential. Results are obtained using calculations with tree level improved gauge and the highly 
improved staggered quark (HISQ) actions with almost physical light and strange quark masses at 
three different values of the lattice cut-off. Our choice of parameters corresponds to a value of 
160 MeV for the lightest pseudo scalar Goldstone mass and a physical value of the kaon mass. 
The three diagonal charge susceptibilities and the correlations among conserved charges have been 
extrapolated to the continuum limit in the temperature interval 150 MeV < T < 250 MeV. We 
compare our results with the hadron resonance gas (HRG) model calculations and find agreement 
with HRG model results only for temperatures T < 150 MeV. We observe significant deviations in 
the temperature range 160 MeV < T < 170 MeV and qualitative differences in the behavior of the 
three conserved charge sectors. At T ~ 160 MeV quadratic net baryon number fluctuations in QCD 
agree with HRG model calculations while, the net electric charge fluctuations in QCD are about 
10% smaller and net strangeness fluctuations are about 20% larger. These findings are relevant to 
the discussion of freeze-out conditions in relativistic heavy ion collisions. 
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I. INTRODUCTION 

The low energy runs currently being performed at the Relativistic Heavy Ion Collider (RHIC) [l| aim at an ex- 
ploration of the QCD phase diagram at non-zero temperature (T) and baryon chemical potential (/is) through the 
measurement of fluctuations of conserved charges, e.g., net baryon number, electric charge and strangeness. For the 
former, first results have been published by the STAR collaboration [2j and preliminary results on the latter two have 
been presented at conferences. A central goal of these experiments is to search for the existence of the QCD critical 
point, a second order phase transition point, that has been postulated to exist at non- vanishing baryon chemical 
potential in the T-fiB phase diagram of QCD [3, [j|. Fluctuations in conserved charges can probe this critical point, 
the cndpoint of a line of first order phase transitions that extends to large baryon chemical potential at non-zero 
quark masses. More generally, the study of fluctuations at any value of the baryon chemical potential probe thermal 
conditions in a medium and provide information on the critical behavior of QCD [5[ . 

At zero baryon chemical potential, already, the analysis of fluctuations of conserved charges and their higher order 
cumulants provides important information about the relation between the QCD chiral phase transition at vanishing 
light quark masses, the cross-over temperature at physical quark masses and the freeze-out conditions observed in 
heavy ion experiments [5j-l7(. Thus, calculations of conserved charge fluctuations at /j,b = will provide unique 
information on freeze-out conditions at the LHC where the baryon chemical potential is small, hb/T sw 0.05. The 
quadratic fluctuations of the net baryon number characterize the width of the probability distribution which has been 
measured at RHIC |2j] and has recently been analyzed in the framework of the hadron resonance gas (HRG) model 8| . 
Although the HRG model provides a rather satisfactory description of global hadron yields at chemical freeze-out 9| , 
its ability to describe detailed properties of strongly interacting matter such as fluctuations of conserved charges and, 
in particular, their higher order cumulants is not obvious. The foundation for the HRG model is given by the Dashen, 
Ma, Bernstein theorem |I2| , which shows that the partition function of strongly interacting matter can be described 
by a gas of free resonances, if the system is sufficiently dilute and the resonance production is the dominant part of 
the interaction among hadrons [13l . |14| | . At very low temperature and high baryon number density, where nonresonant 
nucleon-nucleon interactions become important, as well as at high temperature where strongly interacting matter 
undergoes a transition to the quark-gluon plasma regime and partonic degrees of freedom become dominant, the HRG 
model is expected to be a poor approximation to the thermodynamics of strongly interacting matter. To what extent 
the HRG model provides a good description of strongly interacting matter needs to be explored in detail by comparing 
model calculations with first principal (lattice) QCD calculations. The latter provides the complete description of 
QCD thermodynamics at all values of the temperature and ultimately should set the standard for the interpretation 
of experimental results on strong interaction thermodynamics. Unfortunately, this is, at present, not fully possible. 
For instance, the only satisfactory way to specify the thermal conditions in a heavy ion experiment at the time of 
hadronization is through the comparison of experimental data with HRG model calculations [l5| ■ It has recently been 
pointed out that in the future this may be overcome by comparing experimental data on conserved charge fluctuations 
directly with lattice QCD calculations [5j. Also in order to establish such an approach more firmly it is important to 
understand and quantify to what extent lattice QCD calculations and HRG model calculations agree and in which 
temperature regime the latter provides a reasonable approximation to strong interaction thermodynamics. 

Indeed, some deviations in ratios of higher order cumulants of baryon number fluctuations, calculated within the 
HRG model, from experimental results for cumulants of net proton fluctuations have been observed |£|. Whether 
these deviations can be accounted for within QCD or are of more technical origin related to the restricted phase space 
in which experimental observations have been performed is an open question. In any case, a more detailed analysis 
of the thermal conditions achieved in heavy ion experiments is important. Lattice QCD calculations of fluctuations 
of conserved charges in equilibrium thermodynamics provide a base line for such discussions. Studies of fluctuations 
and higher order cumulants [16[ may reveal differences between HRG model calculations and QCD thermodynamics 
that will appear close to criticality in the QCD phase diagram. 

In this paper we will use the HRG model in its simplest version, i.e. as a sum of non-interacting, point-like particles. 
This is known to provide an accurate description of a dilute, strongly interacting hadron gas [K|]. The need for taking 
into account residual interactions, for instance through the introduction of an intrinsic size of the hadrons [lOj has been 
discussed. It has also been noted that the inclusion of higher mass resonances and an improvement in the strangeness 
sector of the HRG model may be needed to adequately describe pronounced features of hadron production such as 
the enhancement in the K/n ratio [ll(. The advantage of the simplest version, however, is that the HRG model in 
this form is parameter free, while any improvement on this model will introduce further parameters without bringing 
us closer to the actual underlying theory, QCD, in a controlled way. 

Quadratic fluctuations of conserved charges are closely related to quark number susceptibilities [171 ] . Fluctuations of 
net baryon number, electric c harg e and strangeness, as well as correlations among them, have been analyzed in previous 
lattice QCD calculations |18l - [2l| and have also been used to characterize properties of the relevant thermodynamic 
degrees of freedom at low as well as high temperature [161 |22| . |23| . The generic forms of their temperature dependence 



and their scaling properties are understood in terms of universal properties of the QCD partition function and its 
derivatives in the vicinity of the QCD chiral phase transition [7|, l20( . To make use of this knowledge in a quantitative 
comparison with experimental results, lattice QCD calculations close to the continuum are needed. 

In this paper we present an analysis of fluctuations in, and correlations among, conserved charges using numerical 
calculations in (2+l)-flavor QCD at three values of the lattice cut-off Q. For these calculations we exploit an 0{a 2 ) 
improved action consisting of a tree-level improved gauge action combined with the highly improved staggered fermion 
action (HISQ/tree) [2y, [27|. We discuss the cut-off dependence of our results in different temperature intervals 
and consider two different zero-temperature observablcs for the determination of the temperature scale used for 
extrapolations to the continuum limit. This allows us to quantify systematic errors in our calculation. In an appendix, 
we discuss the relation between temperature scales deduced from different zero-temperature observables and the 
propagation of their cut-off dependence into the cut-off dependence of thermodynamic observables. 

II. FLUCTUATIONS OF CONSERVED CHARGES FROM LATTICE QCD; THE HADRON 
RESONANCE GAS AND THE IDEAL GAS LIMIT 

To calculate fluctuations of baryon number (B), electric charge (Q) and strangeness (S) from (lattice) QCD we 
start from the QCD partition function with non-zero light (/i„, fid) and strange quark ()j, s ) chemical potentials. The 
quark chemical potentials can be expressed in terms of chemical potentials for baryon number (hb), strangeness (fJ-s) 
and electric charge (hq), 

1 2 

/"« = -fJ. B + g/^Q , 

1 1 
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Ms = g/^3 - 3^Q - MS ■ (1) 

The starting point of the analysis is the pressure p given by the logarithm of the QCD partition function, 
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Fluctuations of the conserved charges and their correlations in a thermalized medium arc then obtained from its 
derivatives evaluated at (l = {/j.b, Hq, Us) = 0, 
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with fix = Hx/T and X, Y = B, Q, S. Explicit expressions for the calculation of these susceptibilities in terms of 
generalized light and strange quark number susceptibilities arc given in [20| . 

As all these derivatives are evaluated at /2 = 0, the expectation values of all net charge numbers SNx = Nx — N x , 
with Nx (N x ), denoting the number of particles (anti-particles), vanish, i.e., (SNx) = 0. The susceptibilities, i.e., 
the quadratic fluctuations of the charges, arc then given by 

if = ((SN X ) 2 )/VT 3 . (5) 

A. The hadron resonance gas 

We will compare results for fluctuations and correlations defined by Eqs. (|3|) and (Q| with hadron resonance gas 
model calculations. The partition function of the HRG model can be split into mesonic and baryonic contributions, 
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Preliminary results of this work had been presented at Quark Matter 2011 [2J and PANIC 2011 [25 
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TABLE I. Ideal gas values for off-diagonal, Xri = Xii Y ' /T 2 , and diagonal susceptibilities, x 2 — Xu ■ Here X, y = S, Q, 5. 
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where the partition function for mesonic (M) or baryonic (B) particle species i with mass Mi is given by, 
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Here upper signs correspond to mesons and lower signs to baryons; e% = yfk 2 -\- Mf denotes the energy of particle i, 
di is its degeneracy factor and its fugacity is given by 



Zi = exp ((BifiB + QiHQ + SiHs)/T) 



(8) 



With these relations it is straightforward to calculate susceptibilities and charge correlations in the HRG model using 
Eqs. © and (gj. 

We note that a HRG model is defined by specifying the resonance spectrum used to construct the partition function 
in Eq. ©. We use all hadron resonances with masses Mh < 2.5 GeV listed by the particle data group (PDG) in their 
2010 summary tableqj [2i§. It is similar to that used, for example, in Ref. 11 1. One question, that we will discuss in 
the comparison of lattice QCD results with HRG model calculations, is to what extent the strangeness sector is well 
represented in the HRG model calculations. This question has also been addressed in (ll| . 



B. The ideal gas limit 

In the infinite temperature limit, the grand canonical QCD partition function reduces to that of an ideal gas of 
quarks and gluons. In this limit, quark mass effects, including those in the strange quark sector, are negligible and 
we may compare our (2+l)-flavor QCD calculations with a free quark-gluon gas of 3-flavor QCD (Stcfan-Boltzmann 
(SB) gas). This is given by J2<| 
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where the first two terms give the contributions of the gluon and the quark sectors for vanishing chemical potentials. 
After expressing the flavor chemical potentials in terms of /ib , MQ an d Ms as given in Eq. (JTJ , it is straightforward to 
read off the ideal gas values for diagonal and off-diagonal susceptibilities. These are listed in Tabic |TJ 

III. LATTICE CALCULATIONS 

In order to analyze fluctuations of conserved charges, we perform calculations on gauge field configurations generated 
in our study of the finite-temperature transition in (2+l)-flavor QCD [27J |. These calculations were performed with 



2 In this summary table there are a few three starred resonances listed which do not have a known spin assignment. For these we use the 
minimal degeneracy factors. We also checked that the inclusion of some known heavier resonances as well as the inclusion of charmed 
hadrons docs not alter the picture presented here. Moreover, we checked the stability of HRG results by reducing the mass cut-off from 
2.5 GeV to 2.0 GeV. This alters the relevant observables discussed here by at most 2% at T = 200 MeV. 



the HISQ/tree action for three values of the lattice cut-off corresponding to lattices with temporal extent N T — 6, 8 
and 12 and a spatial lattice extent N a = AN T . The calculations cover a temperature range, 130 MeV <T< 350 McV. 
In this temperature range the line of constant physics is defined by tuning the strange quark mass m s to its physical 
value and setting the light quark masses to mi = m s /20, which correspond to a pion mass M^ ~ 160 MeV. 

It is important to note that the definition of physical quark and/or pion masses is not straightforward in calculations 
with staggered fcrmions at non-zero values of the lattice spacing due to taste symmetry breaking. Taste symmetry 
breaking, which is a consequence of the doubling problem in the staggered formulation, gives rise to sixteen pseudo- 
scalar mesons corresponding to the sixteen elements of the Clifford algebra, of which only one, with taste matrix 75, 
behaves as a Goldstone particle at finite lattice spacing. Even in the chiral limit, the other 15 modes receive masses 
of 0(a 2 ) which vanish only in the continuum limit. Thus, as a result of taste symmetry breaking, these sixteen modes 
contribute to observables with different masses. The same is true of all other states, but the problem is most severe 
for Goldstone modes. The size of these effects has been quantified for the staggered fermion discretization scheme 
exploited here (HISQ/tree) by measuring the masses of the sixteen taste pions and by defining a root-mean-squared 
mass [27J |. 

M™ s = \ y/M^ + A/2 o75 + 3A/2 i75 + 3A/£ 7 . + 3Af * 70 + 3 A/2 + M* o + M* . (10) 

Taste symmetry breaking also affects the generation of the ensemble of gauge configurations on which measurements 
are made. To simulate the desired number of flavors, the fourth-root of the staggered fermion determinant is taken 
for each flavor. While this "rooting" trick corrects for the number of flavors in the continuum limit, at finite lattice 
spacings, the masses of states contributing to the partition function are not degenerate. 

With the improved sta gge red fermion action (HISQ/tree) used in our calculation, taste violations, while strongly 
suppressed, are still large [27|. In fact, our analysis showed that the HISQ/tree action has the smallest taste violations 
compared to the other improved staggered actions (p4, asqtad and stout) used in finite temperature calculations. 
Nevertheless, at values of the cut-off corresponding to the transition region (T ~ 160 McV) on the three different 
N T lattices analyzed by us, the RMS masses vary from M^ MS ~ 215 MeV on our finest lattices (N T = 12) to 
MRMS ^ 415 Me y on tne coarses t lattice (7V T = 6), in contrast to the Goldstone pion mass 140 MeV [27J]. We also 
find that, to a good approximation, the difference of M^ M and the physical pion mass M v is proportional to the 
square of the lattice cut-off at a fixed value of the temperature, (aT) 2 ~ 1/N 2 . In this paper we show that the cut-off 
effects resulting from a much heavier RMS mass influence most strongly the electric charge fluctuations as these are 
most sensitive to the contributions from pions. It is worth pointing out, for comparison, that to obtain an RMS mass 
of about 215 MeV achieved with the HISQ/tree action on N T = 12 lattices will require N T s=s 20 with the asqtad and 
stout actions [2lll27j. 

We typically analyzed 10,000-20,000 gauge field configurations per parameter set for 7V r = 6, 10,000-14,000 con- 
figurations for N T = 8 and up to 6000 configurations for N T = 12 lattices. Measurements of all operators needed to 
calculate quadratic fluctuations have been performed every 10 hybrid Monte Carlo time units for N T = 12 and every 
10 or 20 time units for N T = 6,8 lattices. Most of our quark number susceptibilities on the N T = 6 and 8 lattices 
were calculated on GPU-clusters and we used 500-1500 random source vectors for the analysis. The N T = 12 data 
were analyzed using 400 random source vectors at low temperatures and 100 at high temperatures. 

When performing lattice QCD calculations at non-zero temperature, we have to control (at least) two different 
sources of cut-off errors. On the one hand there is the intrinsic cut-off dependence of the observables calculated at 
non-zero temperature at a certain value of the cut-off a -1 . We reduce these by working with tree-level 0{a 2 ) improved 
actions in the gauge as well as the fermion sector. This improvement also insures that at high temperature the cut-off 
dependence of the basic operators entering our calculations is small. An additional cut-off dependence arises due 
to the choice of the zero temperature observable used to set the scale for all finite temperature measurements. We 
investigate two different scale-setting observables: the length scale n extracted from the slope of the static quark 
potential, 

r 2 ^l) =1.0, (11) 

dr /r=n 

and the kaon decay constant fa- The temperature in these units is T ri n = ri/(aN T ) and Tf K /fx = l/(fK<iN T ). 
The efficacy of these observables in setting the scale is complementary in many respects [27| . In particular, we note 
that n has a mild dependence on the quark masses and is well defined even in the infinite quark mass limit. The 
kaon decay constant, on the other hand, has a sizeable quark mass dependence, and at lattice spacings used in this 
study, it is sensitive to taste symmetry violations in the hadron spectrum. To convert to physical units, we use 
n = 0.3106(20) fm Ho| and the latest PDG value for the kaon decay constant, f K = 156.1/\/2 McV. In Fig. [Q wc 
show the difference in the two estimates of temperature over the interval relevant to the calculations performed here. 
The leading order correction contributing to this difference is 0(g 2 a 2 ). 
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FIG. 1. Difference in temperature scales obtained from calculations of n and fx at values of the cut-off relevant for the 
temperature range explored in the calculation of quark number susceptibilities. 
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FIG. 2. Net baryon number fluctuations versus temperature. The left hand figure shows results using the potential shape 
parameter n to set the scale for the temperature. The right hand figure shows the same data using fx to set the scale. Also 
shown are the results obtain from the HRG model and the infinite temperature ideal gas limit (solid lines). 

As noted in J22J , we also find that in the study of charge fluctuations it is advantageous to use the fx scale as it 
absorbs a significant fraction of the cut-off effects, i.e., the cut-off effects are similar and cancel to a large extent in 
the ratio of hadron masses and temperature, M/T. However, it should be stressed that any one observable (r\ or 
fx) cannot eliminate cut-off effects in all observables equally well. We elaborate on this point in more detail in the 
appendix. 



IV. FLUCTUATIONS 



A. Fluctuations in baryon number, strangeness, and electric charge 



We start our discussion of the fluctuations of baryon number, strangeness and electric charge by summarizing 
the data obtained on lattices with temporal extent N T = 6, 8 and 12 in Tables HH IIIII and IIVI and discussing their 
scaling behavior. The continuum extrapolation and comparison to the HRG model as well as the asymptotic high 
temperature ideal gas results are discussed in detail in the next subsection. 

Figure [2] shows results for the baryon number susceptibility with the temperature scale set using t\ (left hand 
panel) and fx (right hand panel). In both cases, we show the results from the HRG model including all resonances 
with mass Mh < 2.5 GeV. The noticeable differences between the left and right hand panels of Fig. [2] are due to 
setting of the temperature scale and we find that the cut-off effects are smaller when a scale based on fx is used. As 
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TABLE II. Quadratic fluctuations of net baryon number, electric charge and strangeness as well as correlations among these 
conserved charges in units of T 2 calculated on lattices with temporal extent N T — 6. We use fx to define the temperature 
scale. 
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-0.0731(18 


0.0252( 7 


0.1081( 6) 


6.445 


165.98 


0.1444(27) 


0.4161(30) 


0.3398(19) 


-0.0919(15 


0.0263( 8 


0.1239(11) 


6.460 


168.57 


0.1610(24) 


0.4287(19) 


0.3623(30) 


-0.1034(17 


0.0288( 5 


0.1294( 9) 


6.488 


173.49 


0.1731(34) 


0.4534(24) 


0.4058(44) 


-0.1170(28 


0.0281( 4 


0.1444(10) 


6.515 


178.36 


0.1948(24) 


0.4794(28) 


0.4606(35) 


-0.1390(19 


0.0279( 5 


0.1608(11) 


6.550 


184.84 


0.2085(17) 


0.5034(19) 


0.5177(30) 


-0.1572(13 


0.0257( 4 


0.1803(10) 


6.575 


189.58 


0.2190(17) 


0.5171(10) 


0.5586(19) 


-0.1719(13 


0.0236( 2 


0.1933( 4) 


6.608 


196.01 


0.2291(23) 


0.5321(18) 


0.6052(39) 


-0.1874(22 


0.0209( 1 


0.2089(10) 


6.664 


207.32 


0.2465(22) 


0.5511(13) 


0.6711(35) 


-0.2123(20 


0.0171( 2 


0.2294( 8) 


6.800 


237.07 


0.2748(14) 


0.5786(14) 


0.7796(27) 


-0.2540(13 


0.0104( 1 


0.2628( 7) 


6.950 


273.88 


0.2901( 7) 


0.5963(10) 


0.8443(18) 


-0.2782( 7; 


0.0060( 1 


0.2831( 5) 


7.150 


330.23 


0.2984( 5) 


0.6073( 9) 


0.8829(13) 


-0.2919( 5; 


0.0032( 1 


0.2955( 4) 



TABLE III. Same as Table \U} but for N T = 8. 



13 T (MeV) 


B /rrn2 

X2/T 


X?/T 2 


S m-,2 
X2/T 


xfi s /r 2 


BQ /rr 2 


X?f/T 2 


6.70 143.25 


0.0404(135) 


0.2741( 64) 


0.1586( 53) 


-0.0227(59) 


0.0088(65) 


0.0679(33) 


6.74 149.03 


0.0764( 67) 


0.3178( 52) 


0.1939( 46) 


-0.0348(40) 


0.0208(32) 


0.0796(25) 


6.77 153.48 


0.1031( 68) 


0.3492( 86) 


0.2365( 59) 


-0.0581(42) 


0.0225(34) 


0.0892(33) 


6.80 158.05 


0.1056( 83) 


0.3721( 46) 


0.2663( 53) 


-0.0632(47) 


0.0212(40) 


0.1015(27) 


6.84 164.31 


0.1520( 42) 


0.4358( 45) 


0.3531( 59) 


-0.0975(30) 


0.0272(20) 


0.1278(22) 


6.88 170.77 


0.1522( 64) 


0.4471( 40) 


0.3808( 66) 


-0.1009(46) 


0.0257(30) 


0.1399(30) 


6.91 175.76 


0.1920( 46) 


0.4943( 36) 


0.4650( 63) 


-0.1361(42) 


0.0280(21) 


0.1645(28) 


6.95 182.59 


0.2034( 55) 


0.5030( 59) 


0.5067( 90) 


-0.1506(55) 


0.0264(28) 


0.1780(46) 


6.99 189.64 


0.2273( 73) 


0.5273( 69) 


0.5779( 87) 


-0.1820(66) 


0.0227(35) 


0.1979(48) 


7.03 196.91 


0.2470( 49) 


0.5483( 36) 


0.6319( 87) 


-0.2006(45) 


0.0232(24) 


0.2157(34) 


7.10 210.21 


0.2547( 33) 


0.5671( 27) 


0.7004(106) 


-0.2213(38) 


0.0167(16) 


0.2395(30) 


7.15 220.15 


0.2650( 25) 


0.5744( 26) 


0.7389( 44) 


-0.2366(32) 


0.0142(12) 


0.2512(24) 


7.28 247.91 


0.2822( 19) 


0.5902( 15) 


0.8063( 37) 


-0.2646(22) 


0.0088( 8) 


0.2708(15) 



TABLE IV. Same as Table M but for N T = 12. 
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FIG. 3. Net strangeness (left) and net electric-charge (right) fluctuations versus temperature obtained from lattice calculations 
at three different values of the lattice cut-off aT = 1/N T . Calculations of fa have been used to fix the temperature scale. Also 
shown are the results obtained from the HRG model and the infinite temperature ideal gas limit (solid lines) as well as three 
short, solid lines for a HRG model in which the pion mass has been replaced by the RMS pion mass relevant for our calculations 
on the N T = 6, 8 and 12 lattices, respectively. 



pointed out above, this feature has also been noted in the analysis of chiral observables [27|, bLfl]. In principle, both 
scales should lead to identical results in the continuum limit, however, the continuum extrapolation is much better 
controlled when the cut-off effects are small. We therefore use the Jk temperature scale in the rest of the paper unless 
stated otherwise. 

In Fig. [3l we show results for fluctuations of net strangeness (left) and net electric charge (right) using Jk to set 
the scale for the temperature. The figure also shows the corresponding HRG results. Both, the strangeness and the 
baryon number (see Fig. [2]) fluctuations agree with HRG results for temperatures below T ~ 160 MeV and show small 
cut-off effects. The electric charge fluctuations deviate significantly from HRG results at all temperatures and show 
large cut-off dependence. These large cut-off effects in xj can mostly be explained as due to discretization errors 
in the lattice hadron spectrum. The dominant contribution to X2 a ^ l° w temperatures comes from pions, while xf 
receives the leading contribution from kaons and %f from nucleons. As discussed in Section IIII1 the pion spectrum 
is strongly affected by taste symmetry violations in the staggered formulation, and it has been shown in earlier work 
that the agreement in jjg between lattice QCD results and HRG model calculations can be improved using a distorted 
spectrum in the HRG analysis [32|, |34j • 

To confirm this observation, we show in Fig. [3] (right panel) the HRG results obtained after replacing the physical 
pion mass by the RMS pion mass relevant to our calculations on the N T = 6, 8 and 12 lattices, respectively. These 
cut-off effects leading to the spectral distortion vary not only with A^- but also with temperatures at fixed N T . We 
therefore parametrize the cut-off dependence of M^ MS , at each N T , using a cubic polynomial fit to the data given 



in [27|]. This allowed us to estimate the modified HRG result as a function of temperature, as shown in Fig. [3]by short 
solid lines (colored lines). The data confirm that below T ~ 155 MeV, the numerical results for jQ are well described 
by this minimally modified HRG model and the major part of the difference is indeed due to the taste symmetry 
breaking effects in the pion sector. Since all other heavier states contributing to the HRG model are assumed to take 
on their physical values, the fluctuations of strangeness and baryon number are not influenced by this modification 
of the HRG model. We discuss these features, together with continuum extrapolations, for the data shown in Fig. [2] 
and [3] in more detail in the following subsection. 

It is worthwhile to clarify our discussion of the comparison of lattice QCD data with the HRG model. When we say 
that the HRG is a good approximation to QCD, we refer to the value of the susceptibility, as is traditional. At the 
temperature where this agreement fails, we observe that not only the value but the slope also deviates significantly. 
Our data also indicate that the curvature starts to deviate 20-30 MeV earlier depending on the observable. These 
derivatives of the susceptibility are related to higher moments, which are increasingly less well captured in the HRG 
analysis and have not been calculated in our lattice simulations. 



B. Continuum extrapolation and approach to the hadron resonance gas estimates 

In this section, we analyze cut-off effects at fixed values of the temperature for three values of the lattice spacing, 
aT = 1/6, 1/8 and 1/12. We perform this analysis for temperature scales defined in terms of both r\ and Jk in order 
to quantify systematic effects at a given lattice spacing and to demonstrate consistency between the estimates in the 
continuum limit. Our simulations on the lattices with the three different N T values have not been done at the same 
values of the temperature. As a result, in order to perform continuum extrapolations at fixed temperature, we have 
used cubic spline interpolations of our data throughout this paper. We propagate errors on the spline parameters. 
In addition to estimating statistical errors based on our entire data sample, we also use the difference of spline 
interpolations performed on two independent sub-sets constructed by choosing every second T-point (even and odd 
points), as an estimate for the systematic errors. The final error on our data is obtained by adding this error estimate 
and the purely statistical error obtained from the full data set in quadrature. 

We first discuss the continuum extrapolation for the three susceptibilities shown in Figs. [2] and [3] at high temper- 
ature, T > 170 MeV. In this regime, cut-off effects are generally small, which, to some extent, is due to the fact 
that our numerical calculations have been performed with an 0(a 2 ) improved action with small cut-off dependence 



of thermodynamic observables in the infinite temperature ideal gas limit |33f . In Fig. [4j we show continuum extrap- 
olations for all three susceptibilities, x 2 > a * t nree representative values of the temperature, T = 170, 190 and 
210 MeV. In each case, and for both temperature scales, r\ and /k, the fits show that the cut-off effects are consistent 
with 0(g 2 (aT) 2 ) corrections and, over the limited range of T, all three susceptibilities can be extrapolated to the 
continuum with an Ansatz that includes corrections linear in 1/N 2 = (aT) 2 . The continuum extrapolated results 
obtained with the two temperature scales agree within errors, and the results obtained on the N T = 12 lattices are a 
good approximation to these. This also is true for lower temperatures. However, in this case extrapolations linear in 
1/N 2 are no longer sufficient for the electric charge and strangeness fluctuations. Systematic effects at 0((aT) 4 ) start 
to become important. This is evident from the data sets at T — 150 MeV, which are shown in Fig. 2] as well. We 
note that also at this temperature, which is the lowest temperature for which we perform continuum extrapolations, 
extrapolations based on the r\ and Jk temperature scales are in good agreement. Having demonstrated consistency 
of the continuum estimate obtained using n and Jk, we, as stated previously, use the scale from Jk in the rest of 
the paper because the slope in the fits is smaller. 

The data for the net charge fluctuations in the temperature interval 120-250 MeV, results of the linear extrapolation 
for xf /T 2 , and quadratic extrapolations for xf/^ 2 an d X2 /T 2 are shown in Figs. [5] and El In Fig. fright) we also 
show the ratio of net baryon number and electric charge fluctuations. The continuum extrapolation shown for this 
quantity has been obtained from the corresponding extrapolations for xf /T 2 and xi /T 2 - 

Continuum extrapolations in the crossover and the low temperature regions require additional considerations be- 
cause the three different conserved charge susceptibilities show different sensitivities to cut-off effects. In order to 
quantify differences from the HRG model results in this temperature regime, and in order to clarify the extent to 
which the HRG model provides a good description of QCD results, we analyze the ratios xi 1X1 ' ■ X — B, Q, and 
S, in Fig. [7] We find that, while baryon number fluctuations start to agree with HRG model results for T<165 MeV, 
the net strangeness fluctuations become larger than the HRG values for temperatures below T ~ 190 MeV and then 
approach the HRG values from above at T<150 MeV. At T ~ 150 MeV, the differences are still (10-20)%. 

The electric charge fluctuations show much larger deviations from the HRG model as is evident from Fig. [7l In 
particular, below T ~ 170 MeV, the cut-off dependence in X2/X2 ^ s l ar g e an d extrapolations including just 

leading order a 2 -corrections fail. As discussed in Section I IV A[ this, to a large extent, is due to the severe cut-off 
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FIG. 4. Continuum extrapolations of net strangeness (xf /T 2 ), net electric charge (xj /T 2 ) and net baryon number (x§ /T 2 ) at 
four values of the temperature using fits linear and, for xf /T 2 and x% jT 2 , at the lowest temperature quadratic in (aT) 2 — 1/N 2 . 
Data at fixed values of the temperature are obtained from cubic spline interpolations. The temperature scale has been 
determined using calculations of ri (circles) and Jk (boxes), respectively. 



dependence of the pion spectrum, i.e., the anomalously large RMS pion mass suppresses fluctuations in the electric 
charge and has a much smaller effect on the baryon and strangeness charges. In short, a continuum extrapolation 
without including the effects of taste symmetry breaking is insufficient. 

The distorted HRG model, which modifies the log of the partition function by replacing M w by M^ MS in the pion 
contribution, exp(— M^/T), however, does describe the data well. In general, the HRG model defined in Eq. [7] suggests 
that, in this temperature regime, cut-off effects in any quantity / may be accounted for by an exponential Ansatz of 
the form 



f(N T ,T) = a(T) + b(T) e - c(T)/N * 



(12) 



which, at high temperatures where cut-off effects become small, reduces to the linear fit in 1/N%, i.e., f(N T ,T) ~ 
a(T) + b(T)/N? and also incorporates the next to leading order quadratic corrections f(N T ,T) ~ a(T) + b(T)/N? + 
c(T)/N*. We, therefore, analyze our data for X2 m * ne transition region, 150 < T < 190 MeV, from the low 
to high temperature phase of (2+l)-flavor QCD using fits linear and quadratic in 1/N% as well as the exponential 
Ansatz given in Eq. (fT2")) . With our current statistical accuracy we are, however, not sensitive to cut-off effects beyond 
0((aT) 4 ). In fact, all our fits performed with the exponential ansatz are consistent within errors with fits based 
on the quadratic ansatz. We, therefore, do not discuss the exponential fits any further in this paper. In the case 
°f X2 we nn d that linear and quadratic fits agree within errors and lead to a x 2 /dof less than unity in the entire 
range of temperatures T > 150 MeV. We therefore use the linear fits to perform continuum extrapolations for xs- 
For strangeness and electric charge susceptibilities we use quadratic extrapolations in the entire temperature range, 
although, as discussed above, we do not observe systematic differences between linear and quadratic extrapolations for 
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FIG. 5. Fluctuations of net strangeness (left) and electric charge (right) in units of T 2 . Calculations of fx have been used to 
fix the temperature scale. Also shown are continuum extrapolated results taking into account cut-off effects up to quadratic 
order in 1/N 2 . The HRG model result and the SB limit is given by the solid lines. 
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FIG. 6. Net baryon number fluctuations in units of T 2 (left) and the ratio of net baryon number and net electric charge 
fluctuations (right). Calculations of fa have been used to fix the temperature scale. Also shown are results from a continuum 
extrapolation taking into account 0(a 2 ) corrections. For \§ IX2 we a ^ so show the ratio of continuum extrapolations constructed 
f° r X2 /T 2 and X2 /T 2 separately. The HRG model result and the SB limit is given by the solid lines. 



T>170 MeV. Using the latter for our continuum extrapolations, however, leads to more conservative error estimates. 

A summary of our continuum extrapolations for X2 IX2 ' an( i the data in the low temperature region is also 

shown in Fig. Due to taste symmetry breaking, the data show significant dependence on N T for T<170 MeV. To 
understand this cut-off effect we compare, in the bottom panel in Fig. [71 results for X2 /^2 with a modified HRG 

model in which the physical pion mass is replaced by (i) the N T dependent RMS mass, HRG(M^ MS )/HRG(AfP h y sical ), 
and (ii) a pion mass of 160 MeV, HRG(Af w = 160 MeV)/HRG(MP hysical ), corresponding to the light quark mass 
actually used in our calculations. Since the line showing HRG^/tt = 160 McV)/HRG(MP hyslcal ) is much closer to 
unity compared to the other three, we confirm that the errors due to simulating at this slightly heavier pion mass are 
significantly smaller than the cut-off effects leading to a much heavier RMS mass. 

In the interesting temperature range relevant to the discussion of freeze-out conditions in heavy ion collisions, 
160 MeV<T<170 MeV, we find that the continuum extrapolated electric charge fluctuations are (10-20)% smaller 
than even the modified HRG model calculation with M w = 160 MeV. For temperatures T<150 MeV, the X2 data 
start to agree with the modified HRG results with M OT = M^ MS and continuum extrapolations using the quadratic 
Ansatz start to agree with the HRG result. 

Strangeness fluctuations on the other hand, both for the N T = 12 data and the continuum extrapolated values, stay 
systematically above the hadron resonance gas result in the temperature range 150-190 MeV. We, therefore, expect 
this feature to survive the continuum extrapolation. Below T<150 MeV, fluctuations in the strangeness charge show 
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FIG. 7. Deviations of X2' from HRG model results. Data on N T — 6, 8 and 12 lattices are shown using triangles, boxes 
and circles, respectively. Thin lines are the result of a cubic spine interpolation of the data at fixed N T . The broad bands 
are the result of the continuum extrapolations. In the bottom panel, the thick black line just below unity is the ratio 
HRG(M,r = 160 MeV)/HRG(M£ hysical ). The other thick lines show the ratio HRG(M* MS )/HRG(M£ hysical ) for each of the 
three N T lattices. 



an N T dependence, which is most likely again due to taste symmetry violations. More data are required to study this 
issue further. 

Net baryon number fluctuations arc consistent with HRG model results for temperatures below T ~ 160 MeV, 
although statistical errors on our N T = 12 data set make the quantification of possible deviations from the HRG 
result in this temperature range difficult. For larger values of the temperature, the estimates and errors in them grow 
progressively smaller. Our linear extrapolations suggest that in the temperature interval 160 MeV<T<170 MeV, 
deviations from the HRG model calculations are at most 10%. A confirmation of this in our data, through the 
inclusion of quadratic corrections, however, requires better statistics. Data obtained with the stout action [2l| also 
suggest that the xf stays close to the HRG model result in this temperature range. 

In order to reflect the influence of systematic effects on our continuum extrapolation we varied fit ranges and 
distribution of knots in the smooth spline interpolation. Moreover, in order to account for possible underestimates of 
errors on the individual data points, we divided our data samples into two independent sub-sets consisting of even and 
odd T-values by rank-order. We used the differences in these fits and the fit to the full data sample as an additional 
error on our spline interpolations. 

To summarize, our continuum extrapolated values of X2 l^ 2 are given in Table [Vj Extrapolations of X2 l^ 2 
used the Ansatz, Eq. (fl2)) truncated at 0(1/N*), i.e., we used the quadratic extrapolations, whereas for xf /T 2 the 
exponential was truncated at 0(1/N%), i.e., we used the linear extrapolations. For X2 l^ 2 we also find that the 
quadratic extrapolated values agree with the exponential Ansatz. Our extrapolated results are in good agreement 
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TABLE V. Continuum extrapolated results for the quadratic fluctuations of net baryon number, electric charge and strangeness 
densities and the correlations among them. Results for X2 /T and \2 /T are obtained from quadratic fits and those for \2 /T 
from linear fits. 

with the recently published analysis using the stout action [21|. 

Lastly, in Fig. [5] (right), we show the ratio of the net baryon number and electric charge fluctuations, X2 1x2 ■ ^ 
approaches the HRG model result from above and starts to agree with it for T<150 MeV. The continuum extrapolation 
here is based on the linear and exponential extrapolations for xf an d X2 > respectively. In the region of interest to 
heavy ion phenomenology, this ratio varies between 0.29(4) (at 160 MeV) and 0.35(4) (170 MeV). Thus, fluctuations 
in net electric charge could be 3-4 times larger than in the net baryon number in the vicinity of the freeze-out 
temperature. 



CORRELATIONS 



Probes of the structure of QCD at finite temperature include correlations among different conserved charges. These 
correlations show characteristic changes in the crossover region between the low and high temperature phases of QCD, 
which are correlated with changes in the relevant degrees of freedom. They also provide insight into the applicability 
of HRG model calculations at low temperatures. The change in correlations between baryon number and electric 
charge, Xu > 1S expected to be particularly striking as one goes from the low to the high temperature phase. At low 
temperatures this correlation is dominated by the contribution of protons plus anti-protons. Consequently, within 
the HRG model it rises exponentially with temperature in this region. In the high temperature limit of (2+l)-flavor 
QCD, however, Xii vanishes as the quarks become effectively massless, (toj/T — > 0), and the weighted sum of the 
charges of up, down and strange quarks vanishes. Results for X11 shown in Fig. [8] arc consistent with this picture. 

The correlations of strangeness with the baryon number and the electric charge, Xii S an d X11 j are sensitive to 
changes in the strangeness degrees of freedom [22], [33, [36| . Results for the temperature dependence of these correlations 
are shown in Fig. [5J They approach the Stefan-Boltzmann value, 1/3, of a massless three flavor quark gas at high 
temperatures. As observed in the case of the quadratic fluctuations, on decreasing the temperature towards the 
transition region, these correlations first overshoot the HRG model result and then approach HRG value from above 
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FIG. 8. Correlations of electric charge with baryon number versus temperature. The temperature scale has been set using fa- 
The solid line shows the result for the HRG model. The band shows continuum extrapolations that take into account cut-off 
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FIG. 9. Correlations of strangeness with baryon number (left) and electric charge (right). In both cases, fx has been used to 
set the scale. The solid lines show the result for the HRG model. For \i\ /T 2 the band shows continuum extrapolations that 
take into account cut-off effects linear in 1/N 2 while for Xn /T 2 also quadratic corrections have been accounted for. 



at about 150 MeV. This overshoot is more pronounced for —Xii than for X\\ . Also shown in Fig. [9] are continuum 
extrapolations (bands) which in the baryon sector, i.e., for B-S and B-Q correlations, are based on fits linear in 1/N%, 
whereas in the meson sector, i.e., for Q-S correlations, quadratic corrections are also taken into account. This reflects 
the larger sensitivity of the latter to taste violations that also has been observed for the quadratic strangeness and 
electric charge fluctuations. 

In the isospin symmetric case considered in this study, the flavor correlations Xii and Xn are equal. Also, the two 
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QS 



correlations 2xn and X\i arc related to each other through the quadratic strangeness fluctuations, i.e., 2x{{ 



\f 



All 



1361 ] . One can then write the following relationships between the charge correlations and quark- flavor fluctuations: 
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(13) 



At high temperatures Xii ~ — Xii S because Xn receives perturbative contributions only at 0(g 6 ln(l/g 2 )) and is 
therefore small [23.I37J. On the other hand, corrections to x s 2 from the ideal gas limit are dominant as they are 0(g 2 ). 
All three charge correlations show significant deviations from the ideal gas limit even at twice the transition 
temperature (see data in Tables [TTl IIIII and IIV|) . These deviations are due to large contributions of flavor fluctuations, 
such as to xl discussed above. The leading order perturbative correction can be eliminated by forming suitable ratios 
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FIG. 10. Correlations between baryon number and strangeness (left) as well as electric charge (right) normalized to strangeness 
and net baryon number fluctuations, respectively. The solid line shows the HRG model result. 



that can be used to analyze experimental data on charge fluctuations [22| . l35l . \3t 
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At high temperature the deviations from the ideal gas value of unity are now due to Xif/xii f° r example, Cbs = 
1 + 2xii/x2- Data for Cbs is shown in Fig. fTOlf left) and, in comparison to the quadratic strangeness fluctuations 
(Fig. [3]) or baryon number strangeness correlations (Fig. [5]), the approach to the ideal gas limit is much more rapid. 
This shows that the flavor correlation 2x"f/x2 is already small for T>1.2T C . It is, however, large in the vicinity of 
the transition temperature. 

The behavior of the third ratio that one can analyze, Cbq = Xii /X2 ■> ls somewhat different as leading order per- 
turbative corrections do not cancel completely due to differences in the light and strange quark masses. Consequently, 
the approach to the ideal gas limit is slower as can be seen in Fig. [TOT right). 



VI. DISCUSSION AND CONCLUSION 



We have analyzed quadratic fluctuations and correlations among conserved charges in (2+l)-flavor QCD. We find 
that as the temperature is decreased from the high temperature phase, the net baryon number fluctuations start 
to agree with the hadron resonance gas model below 165 MeV while for electric charge fluctuations this happens 
only below 150 MeV. The fluctuations of net strangeness overshoot the HRG model values at T ~ 190 MeV. In the 
temperature range relevant to the discussion of chemical freeze-out in heavy ion collisions, 160-170 MeV, strangeness 
fluctuations are systematically larger than the HRG model result by about 20%. These detailed differences between 
QCD calculations and HRG model results should become manifest when experimental data on the probability distri- 
butions of net charge fluctuations |2j is analyzed. In fact, quadratic fluctuations characterize the bulk structure of 
these distributions which, in the Gaussian approximation at fis = 0, is given by 



P(N X ) = c-^/W 7 ^) , X = B, S, Q 



(15) 



where Tf denotes the temperature at the time of chemical freeze-out and Vf is the freeze-out volume. The results 
presented here suggest that the largest deviations from HRG model calculations occur in the probability distributions 
for electric charge and strangeness fluctuations. We give a summary of the fluctuations in conserved charges for 
temperatures in the transition region in Table IVII 

In the phenomenologically interesting temperature regime 160 MeV<T<170 MeV, continuum extrapolated results 
for %2 ) are smaller than the HRG model results by about 10-20%. Even though at temperatures below 150 MeV, 
estimates for electric charge fluctuations have large systematic errors due to the distortion of the light meson spec- 
trum in all staggered formulations, our analysis shows that these effects can be taken into account when performing 
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T [MeV] 



B i B.HRG v / 
X2 IX2 X2 IX 



Q.HRG S 1 S.HRG 
X2IX1 



BS 
X2 



I BS.HRG 

1X2 



A 2 



? /xf°- HflG 



QS 1 C, 

Xi IX2 



155 
160 
165 
170 



1.049(79) 0.924(36) 1.240(116) 

1.020(80) 0.895(41) 1.235(111) 

0.972(72) 0.861(44) 1.212(106) 

0.898(60) 0.818(41) 1.171( 96) 



1.353(159) 
1.384(135) 
1.356(116) 
1.280( 96) 



0.804(86) 
0.717(67) 
0.633(51) 
0.544(42) 



1.139(74) 
1.144(68) 
1.150(67) 
1.144(64) 



TABLE VI. Quadratic fluctuations of net baryon number (SNb), electric charge (SNq) and strangeness (5Ns) densities and 
correlations among these conserved net charges in the crossover region from the low to high temperature regime of QCD. We 
give results for quadratic fluctuations and correlations calculated in QCD relative to the corresponding HRG model results. 

continuum extrapolations. Wc show that the resulting continuum estimates in the transition region lie below the HRG 
estimates even after corrections accounting for the distorted pion spectrum have been applied. Thus, our conclusion 
is that at the highest RHIC energies and at the LHC, the width of probability distribution for the net electric charge 
should be narrower than what HRG model calculations would suggest since %% in Eq. [15] is smaller. 

In the case of the net baryon number fluctuations, deviations from HRG model results start to become statistically 
significant only for T>165 MeV and will therefore be hard to quantify. We also find that the ratio of the net baryon 
number and electric charge fluctuations, presented in Fig. [BJ approaches the HRG model result from above and starts 
to agree with it for T<150 MeV. In the transition region this ratio is 



ig-~ (0.29- 

X2 



0.35) for 160 MeV < T < 170 MeV 



(16) 



i.e., fluctuations in net electric charge are expected to be about three to four times larger in the vicinity of the freeze- 
out temperature in heavy ion collisions than net baryon number fluctuations. It is worth noting that a comparison of 
X2 7x2 with the ratio of proton to net charge fluctuations, which is accessible in heavy ion collisions, will allow us to 
relate fluctuations in the proton number to the fluctuations of the conserved net baryon number [38[ . 

Finally, we point out that the continuum extrapolated results presented here, and summarized in Table IVII in a 
temperature regime relevant to the freeze-out conditions in heavy ion collisions, do not rely on any uncertainties in 
the determination of the QCD crossover temperature or its characterization through different obscrvablcs as discussed 
in [27j. Systematic errors on the temperature values listed in the first row of Tablc lVII can only come from uncertainties 
in the zero-temperature observable used to determine the temperature scale. We estimate these uncertainties to be 
less than 2 MeV in our calculation [27j |. 
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APPENDIX: CHOICE OF TEMPERATURE SCALE AND THE HADRON SPECTRUM 



In this appendix we discuss the effect of setting the scale using different observables with mass-dimension one 
calculated in zero-temperature simulations. We denote lattice observables measured in units of the lattice spacing 



and calculated at a given j3 by Oi{j3). In the limit f} = 10/.g 2 — > 00, these approach their physical value Of v as 



Oitf) 



o 



phy 



A L 



■R(f3)(l+g%R 2 (P)) 



(17) 
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where Al is the QCD scale, R(/3) = a(/3)Ai is the /3-function, and only the leading correction 0(g 2 a 2 ), relevant to 
our tree-level improved staggered formulation, has been retained. 

Consider using the observable Oi to define the temperature scale Ti. On a lattice of temporal size N T , this 
temperature is given in terms of 0\ v as 

Qphy 

™ " oMK |18 > 

The ratio of any two such temperature scales is then given by 

TiQg) of v o 2 (P) 
T 2 (J3) Oi(/3) ofy 

l+b 2 g 2 R 2 ((3) 
~ l + b l9 2 R 2 (/3) 
~l + (b 2 -b 1 )g 2 R 2 (l3) . (19) 

Using Eqs. (|18|) and (jT^J) we can express the observable 02(/3) in terms of the temperature scale Xi obtained from 0\ 
as 

ftPhy fjphy 

2 ((3)Nr = ^ « ^_ (1 + fe - b l )g 2 R 2 (fi)) . (20) 

J2 J-l 

This shows that if an observable of interest (0 2 ) has a cut-off dependence similar to that of observable (Oi) used to 
determine the temperature scale T\, i.e., b 2 — b\, then 2 ((3)N r , as an estimate of Oj V /T\, has small cutoff effects. 
Since all the bi — b\ need not be small, improving the scaling behavior of one observable does not, in general, imply 
improvement in all observables. 

In the low temperature regime of QCD, the relevant degrees of freedom are hadrons with masses Mh- If the hadron 
resonance gas is a good approximation in this regime, continuum extrapolations of lattice data are better controlled 
if a temperature scale is chosen such that all the lattice estimates of Mjj/T have small cut-off dependences. We find 
that using a hadronic observable such as Jk improves the scaling behavior of the susceptibilities and correlations 
between charges as shown in Fig. [4] 

The pion sector is, however, different and has enhanced cut-off effects due to taste symmetry breaking. One does 
not, therefore, expect to absorb all these effects with a choice of the temperature scale. For this reason we had to 
modify the HRG analysis to compare with lattice data for charge fluctuations which are dominated by contributions 
from the pions. 
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